Source signature determination and multiple reflection reduction

ABSTRACT

The signature of an energy source is determined from an inverse scattering Born series representing multiple reflected energy. The series comprises a polynomial in the inverse of the signature and has recorded data as the first term. The series is truncated, preferably to the first two terms to permit an analytical determination of the signature to be found. The value of the inverse signature which minimises the energy represented by the sub-series is found and this corresponds to the desired source signature. An iterative scheme may be adopted to improve the match to the actual signature so as to take into account the errors caused by truncating the scattering series.

The present invention relates to a method of determining a signature ofa source of energy and to a method of reducing the effects of multiplereflected energy. Such a method may be used in medical imaging andnon-destructive evaluation, where a sample is irradiated with energy,for instance in the form of x-rays or ultrasound, in order to determinethe internal structure of the sample non-invasively. Such a method mayalso be used with marine seismic reflection data obtained by means ofseismic sources and receivers located in water and towed behind aseismic exploration vessel. These methods may be used during actualsurveying and/or subsequently with recorded data from such surveys.

During seismic surveying, a seismic source is repeatedly actuated andseismic receivers, such as hydrophones in marine seismic surveying,receive energy direct from the sources and reflected from variousboundaries or interfaces. In the case of marine seismic surveying,energy propagates into the earth and is reflected back to thehydrophones from subterranean boundaries or interfaces, for instancebetween strata of different types. The recorded seismic data can beprocessed to reveal information about the structure of the earth in thearea being surveyed. However, such reflections are contaminated by otherreflection paths. For instance, energy from the sources is reflected atthe sea surface directly to the hydrophones. Also, energy can bereflected more than once between the sources and the receivers. Suchmultiple reflections can take place within the earth. Also, energyinitially travelling downwards from the sources can be reflectedupwardly and then downwardly again from the surface of the sea beforearriving at the hydrophones. Reflections of this type are referred to as"free-surface multiple reflections". Free-surface multiple reflectionscan be classified according to their order, which is equal to the numberof reflections from the free-surface. Thus, first order free-surfacereflections comprise energy initially travelling downwardly from thesources (as opposed to "ghosting" where energy travels upwardly and isreflected from the free-surface), is reflected upwardly from the sea bedor a boundary below the sea bed, and is then reflected downwardly fromthe free-surface to the hydrophones. Second order free-surface multiplereflections undergo two downward reflections from the sea-surface beforebeing detected by the hydrophones, and so on.

In order to remove or reduce the effects of multiple reflectionsincluding free-surface multiple reflections in seismic data, a goodknowledge of the seismic source signature is required. The seismicsource signature is the wave shape e.g. the pressure waveform withrespect to time, of the seismic energy emitted by the seismic sources.Removal or reduction of the effects of multiple reflections in seismicdata is performed before "stacking", the known process of summingseismic data related to each sub-region of the region being explored soas to improve the effective signal-to-noise ratio of the data. Goodknowledge of the seismic source signature also enhances other dataprocessing techniques such as deconvolution, migration and inversion forelastic parameters. Thus, good knowledge of the seismic source signaturecan significantly enhance the final data processing products.

There are known techniques for determining the signature of a seismicsource. One such known technique relies on a "convolution model" forseismic data, in which the model is defined as the convolution of thesource signature with the impulse response of the earth (including theeffects of reflections, refractions, multiples, and diffractions). Thistechnique requires known information about some part of the sub-surface,such as the sea bottom, so as to allow estimation of the sourcesignature to be treated as a linear inverse problem. However, such sub-surface information is often not available and the technique must thenrely on statistical methods which require more severe assumptions, forinstance about the phase of the signature or the "whiteness" of thereflectivity sequence. Accordingly, the practical use of this techniqueand its accuracy are limited.

Other known techniques comprise actually measuring the signature. Forinstance, in very deep water, the direct wave from the sources to thehydrophone can be measured without interference from reflections fromthe sea bed and sub-strata. "Ghosting" reflections from the sea surfacecan be subtracted with a sufficient degree of accuracy to allow thesource signature to be directly determined. However, it is generally notconvenient to perform such measurements as this generally requirestowing the sources and hydrophones to regions of deep water and thenreturning to a region to be surveyed. Further, the source signature mayvary with time. For instance, where the or each source comprises anarray of individual sources, such as air guns or water guns, one or moreof the individual sources may cease operating during actual surveying,in which case the signature of the whole source will change and willtherefore not be known with sufficient accuracy.

Other techniques for measuring the source signature while performingseismic surveys are also known. However, many of these techniquesrequire special data acquisition geometries i.e. special arrangements ofsources and/or hydrophones, which are often not available. Examples ofsuch techniques are disclosed in EP-A-0 066 423, in Landr, M. andSollie, R., 1992, Source signature determination by inversion:Geophysics, 57, 1633-1640, and in Weglein, A. B. and Secrest, B. G.,1990, Wavelet estimation for a multidimensional acoustic or elasticearth: Geophysics, 55, 902-913.

According to a first aspect of the invention, there is provided a methodof determining a signature of an energy source as defined in theappended claim 1.

According to a second aspect of the invention, there is provided amethod of reducing the effects of multiple reflected energy as definedin the appended claim 11.

Preferred embodiments of the invention are defined in the other appendedclaims.

It is thus possible to provide a technique which may be applied tomarine seismic reflection data to allow source estimation for "prestack"data processing. This technique allows the source signature whichpermits reduction or removal of events due to first order free-surfacereflections to be found. This technique makes use of the formulation ofthe relationship between free-surface reflections and the sourcesignature as a scattering Born series as disclosed in Carvalho, P. M.,Weglein, A. B. and Stolt, R. H., 1991, Examples of a nonlinear inversionmethod based on the T matrix of scattering theory: Application tomultiple suppression: Mtg. Soc. Expl. Geophys., Expanded Abstracts,1319-1322. The series is constructed exclusively with seismic data andthe source signature without any knowledge of the subsurface other thanthe velocity of sea water. Thus, by using the scattering Born series forremoving free-surface multiples, the source signature may be determinedwithout any of the assumptions usually associated with known techniquesbased on the classical convolution model. Further, special acquisitiongeometries of sources and receivers are not required. The sourcesignature may be determined directly from data recorded during seismicsurveying. Further, the source signature thus determined may be used toreduce or remove the effects of free-surface multiple reflections so asto provide a substantial improvement in the quality of the data beforestacking and when other processing techniques are applied.

By using only the first two terms of the inverse scattering series,which is physically equivalent to dealing only with first orderfree-surface reflections, the relationship between free-surfacereflections and the source signature has effectively been made linear asopposed to the nonlinear relationship which applies when higher orderfree-surface reflections and higher order terms of the inversescattering series are considered. The mathematical problem of findingthe source signature which minimises the "energy" of the seismic data(corresponding to substantially removing the effects of first orderfree-surface reflections) has an analytical solution and, in many cases,this solution permits the source signature to be determined to a desireddegree of precision. However, this precision is limited by truncatingthe inverse scattering series and, if desired, a better approximation tothe source signature may be obtained iteratively, for instance byfinding a series of corrections to the source signature whichprogressively reduce the energy of the seismic data. Such iterationsthus correct at least partially for the effects of truncation of theinverse scattering series. This iterative process may be stopped when asufficient degree of precision has been achieved, for instance whenfurther reductions in energy are not considered significant or whenfurther corrections to the source signature are not consideredsignificant.

The invention will be further described, by way of example, withreference to the accompanying drawings, in which:

FIG. 1 is a schematic diagram illustrating marine seismic surveying;

FIG. 2 is a graph illustrating synthetic seismic reflection data as timeagainst offset showing the effect of free-surface multiple reflections;

FIG. 3 is a graph similar to that of FIG. 2 illustrating removal of thefirst order free-surface multiple reflections;

FIG. 4 is a graph similar to that of FIG. 2 illustrating removal offirst and second order free-surface multiple reflections;

FIG. 5 is a graph similar to that of FIG. 2 illustrating removal offirst, second, and third order free-surface multiple reflections;

FIG. 6 is a graph of amplitude against time showing the amplitudespectrum of the actual source signature and the initial estimated sourcesignature corresponding to the data shown in FIG. 2;

FIG. 7 is a graph of amplitude against frequency showing the frequencyspectrum of the actual source signature and the initial estimated sourcesignature corresponding to the data shown in FIG. 2;

FIG. 8 is a graph of phase angle against frequency showing the phasespectrum of the actual source signature and the initial estimated sourcesignature corresponding to the data shown in FIG. 2;

FIG. 9 is a graph similar to that shown in FIG. 2 illustrating theseismic data after removal of free-surface multiple reflections usingthe initial estimated source signature shown in FIG. 6; and

FIGS. 10 to 13 correspond to FIGS. 6 to 9, respectively, and illustratethe estimated source signature and the effects of multiple removal afterfive iterations to reduce the energy of the seismic data.

FIG. 14 is a flow-chart illustrating an example in accordance with theinvention.

For the purpose of illustration, the use of the present invention withmarine seismic reflection data will be described in detail. However,these techniques may be applied to situations where multiple reflectedenergy can be represented by an inverse scattering series.

FIG. 1 illustrates a typical arrangement for performing marine seismicexploration. The drawing shows a section through the earth 1 below thesea 2 with the sea bed shown at 3 and the sea surface shown at 4. Thesea surface 4 constitutes an interface between the air and the water andthus constitutes a free-surface in terms of seismic reflection data. Anexploration vessel 5 tows a seismic source 6, for instance comprising anarray of air guns. The vessel 5 also tows a seismic "streamer" 7, whichis cable several hundreds or thousands of meters in length carryinghydrophones regularly spaced along the cable and connected to datarecording apparatus on board the vessel 5. The spacing between adjacenthydrophones is typically of the order of one or a few tens of meters.

During exploration, the vessel 5 tows the source 6 and streamer 7 atrelatively low speed along parallel lines above the sea bed 3 whileactuating the source 6 at regular intervals. The seismic signalsarriving at the hydrophones of the streamer 7 are recorded on board thevessel with or without data processing. The recorded seismic data aregenerally subjected to further processing elsewhere so as to revealinformation about the structure of the earth 1 below the explorationregion.

The energy emitted by the source 6 arrives at the hydrophones of thestreamer 7 via various different types of propagation paths. Thepropagation paths which are required for processing of the seismic datato reveal information about the structure of the earth are downwardpaths such as 8 followed by upward paths such as 9 to a typicalhydrophone 10 following a single reflection, for instance at a boundarybetween different substrata illustrated at 11. However, these reflectiondata are contaminated with data arriving via other propagation paths.For instance, there is a direct path from the source 6 to eachhydrophone, for instance shown as the path 12 to the hydrophone 13.Also, some energy initially travels upwardly from the source 6, forinstance shown by the path 14, to be reflected downwardly by the seasurface 4 as shown by the path 15. The energy reflected downwardly bythe sea surface 4 also travels to the sea bed 3 and into the earth 1 tobe reflected back to the hydrophones. This gives rise to a virtualsource shown in broken lines at 16. This effect is known as "ghosting".

In general, direct waves following propagation paths such as 12 andvirtual sources 16 resulting from reflection at the sea surface 4 do notcause a problem as their effects can relatively easily be subtractedfrom the seismic signals received by the hydrophones, such as 10 and 13.

The seismic signals are also contaminated by multiple reflections andother effects such as refractions and diffractions. Multiple reflectionscan occur within the earth 1. However, the present invention isconcerned with free-surface multiple reflections i.e. energy from thesource 6 which initially travels downwardly but is received by thehydrophones after at least one downward reflection from the sea surface4. The propagation of a typical free-surface multiple reflection isshown in FIG. 1. Energy from the source 6 travels initially downwardlyalong a path 17. Part of the energy is reflected by the sea bed 3 tofollow the path 18. This energy is substantially wholly reflected by thesea surface 4 along a path 19 and then again partially reflectedupwardly by the sea bed 3 along a path 20 to arrive at the hydrophone10. Such free-surface multiple reflections are categorised by the numberof downward reflections which take place at the free-surface or seasurface 4 following the initial downward propagation of the energy fromthe source 6. Thus, a first order free-surface multiple reflectionundergoes one downward reflection at the surface 4 as shown in FIG. 1. Asecond order free-surface multiple reflection undergoes two downwardreflections at the surface 4, and so on. The present invention providesa technique for determining a signature of the source 6 and provides atechnique for removing or reducing the effects of free-surfacemultiples, as will be described hereinafter.

In order to illustrate this technique, the derivation of the method willbe described with respect to a theoretical two-dimensional earth modelcomprising an inhomogeneous solid half-space overlain by a homogeneousfluid (water) layer. The sea surface 4 constitutes a free-surfacebecause acoustic pressure vanishes in order to avoid infiniteacceleration of the surface layer. The reflection coefficient at thisinterface is almost equal to -1. It is assumed that the signal recordedfrom a typical hydrophone has had the effects of the direct wave andghosting removed. Referring to horizontal (x) and vertical (z) spatialcoordinates and a time coordinate t, the source position is denoted by(x_(s),z_(s)) and the hydrophone (receiver) position is denoted by(x_(g),z_(g)). The seismic data representing the pressure variation at(x_(g),z_(g)) and at time t for a source located at (x_(s),z_(s)) isgiven by D₀ (x_(s),z_(s) t; x_(s),z_(g)). By resetting the time variableto zero for each source and simplifying such that z_(g) =z_(s) =0, theseismic data may be equivalently represented by D₀ (x_(s),x_(g),t). Thisis Fourier-transformed to the (w-k) domain with respect to x_(s),x_(g),tto give D₀ (k_(s),k_(g),w), where k_(s),k_(g), and w are theFourier-transformed variables corresponding to x_(s),x_(g), and t,respectively. To simplify notation, the Fourier-transformed data will berepresented hereinafter by D₀, and corresponding simplified notationwill be used throughout.

As disclosed by Carvalho, P. M., Weglein, A. B. and Stolt, R. H., 1991,Examples of a nonlinear inversion method based on the T matrix ofscattering theory: Application to multiple suppression: Mtg. Soc. Expl.Geophys., Expanded Abstracts, 1319-1322, the scattering Born series forremoving events due to free-surface reflections i.e. free-surfacemultiple reflections, from the seismic data D₀ can be written as:

    D.sub.p =D.sub.0 +AD.sub.1 +A.sup.2 D.sub.2 +              (1)

where D_(p) is the data without free-surface multiples and A=A(w) is theinverse of the Fourier-transformed source signature S(w), which isassumed to be only time dependent and not to vary with angle or sourceposition i.e.

    A=1/S.

The pressure fields D₁,D₂, . . . are given by ##EQU1## where n is apositive integer and ##EQU2##

The constant c is the seismic velocity of water and k is a generichorizontal wave number. The Born series in equation (1) for removingfree-surface multiples is as follows. The series is constructed usingthe (Fourier-transformed) seismic data D₀ and the inverse source A only.The first term D₀ of the series is the actual data. The second termcontaining D₁ removes first order free-surface multiples. The next termincluding D₂ removes second order free-surface multiples, and so on.Thus, in order to remove the effects of free-surface multiplereflections, knowledge of the source signature, and hence its inverse A,is required.

In order to determine the source signature (via its inverse A), it isnecessary to find a solution to the equation (1), which is a polynomialin A. Equation (1) is a polynomial of high order and finding itssolutions is therefore a difficult nonlinear problem. However,truncating of the polynomial to its first two terms corresponds toremoval of first order free-surface multiples and yields a polynomial offirst order which can be solved analytically. Thus, truncating equation(1) to its first two terms gives:

    D.sub.f =D.sub.0 +AD.sub.1 +ε.sub.T                (4)

where D_(f) are the data from which the first order multiples have beenremoved and e_(T) describes the effects due to truncation, is small, andis nonlinearly related to the inverse source A.

Removing the effects of first order free-surface multiples reduces theacoustic energy represented by the seismic data. When the sourcesignature is correctly determined, the reduction in energy will bemaximized. Thus, it is required to find the source signature for whichthe energy of D_(f), i.e. the data after the removal of first ordermultiples, is minimum. By using the known least squares technique, thesignature can be defined by that A which minimizes:

    S(A)=∥D.sub.f ∥.sup.2 +∥A∥.sup.2(5)

where:

    ∥D.sub.f ∥.sub.2 =∫dk.sub.g ∫dk.sub.s ∫dω.D.sub.f.W.sub.D.D*.sub.f                   (6)

and

    ∥A∥.sup.2 =σ.sup.2 ∫dw∫dw'.A(ω).W.sub.A.sup.-1 (ω,ω').A*(ω')(7)

The asterisk denotes complex conjugate. W_(D) =W_(D) (k_(s),k_(g),w) isa weighting function describing errors in the data and W_(A) (w,w')describes a priori information about the source. The term ∥A∥² isintroduced to guarantee the stability of the solution. To simplifysubsequent inversion formulae, the constant s² has been introduced intothe definition of ∥A∥². This criterion corresponds to finding the sourcewhich minimizes the energy of data after the removal of the first ordermultiples. The source signature which permits the removal of the firstorder multiples will also permit the removal of higher order multiplesso that, by finding the source signature in accordance with thistechnique, it may be used in equation (1) to remove free-surfacemultiples up to any desired order.

The a priori information described by the weighting function W_(A)(w,w') may be derived, for instance, from a deterministic estimate ofthe source signature, for instance using the technique disclosed inEP-A-0 066 423. Thus, a priori knowledge such as smoothness of thesource spectrum may be incorporated through correlation betweenfrequencies.

In order to solve equation (5) so as to find the (inverse) sourcesignature A, an initial analytical solution to equation (5) can be found(by ignoring the truncation errors (ε_(T))) by inversion to give:##EQU3## where

    N(ω)=∫dk.sub.g ∫dk.sub.s.D.sub.0.W.sub.D.D.sub.1 *(9)

and

    Q(ω)=∫dk.sub.g ∫dk.sub.s D.sub.1.W.sub.D.D.sub.1 *(10)

and A.sup.(0) is the initial solution.

In some circumstances, the initial solution will comprise a sufficientlyaccurate and precise determination of the source signature so that itmay then be substituted into equation (1) so as to provide data D_(P) inwhich free-surface multiple reflections have been attenuated to asufficient degree. Also, the signature thus determined may be used forother purposes, for instance during subsequent data processing. However,where a more accurate determination of the source signature is required,this may be achieved by setting up an iterative scheme, of whichA.sup.(0) is the starting solution, in order to accommodate thetruncation errors which were ignored in the analytical solutionrepresented by equation (8). Effectively, the nonlinear part of theproblem is solved through the iterative linear solutions.

Initialisation of the iterative scheme is as follows:

    D.sub.f.sup.(0) =D.sub.0 +A.sup.(0).D.sub.1                (11)

and

    D.sub.0.sup.(1) =D.sub.f.sup.(0)                           (12)

In general, A.sup.(0) permits a significant reduction of first ordermultiple energy through D_(f).sup.(0). The iterations are set up to findthe correction to A.sup.(0) which allows the removal of the remainingfirst order free-surface multiple energy.

For each D₀.sup.(n) representing data containing the residual energy ofthe first order multiples, D₁.sup.(n) can be found from equation (2).The correction δA.sup.(n) to A.sup.(n-1) can then be found by minimizing

    S(δA.sup.(n))=∥D.sub.f.sup.(n) ∥.sup.2 +∥δA.sup.(n) ∥.sup.2              (13)

This equation can be solved analytically in accordance with equation(8). D₀ is then replaced by D₀.sup.(n) and D₁ is replaced by D₁.sup.(n).Finally, the source signature is updated in accordance with:

    A.sup.(n) =A.sup.(n-1) +δA.sup.(n)                   (14)

and the data are updated in accordance with

    D.sub.f.sup.(n) =D.sub.0.sup.(n) +δA.sup.(n).D.sub.1.sup.(n)(15)

    D.sub.0.sup.(n-1) =D.sub.f.sup.(n)                         (16)

This iterative process may be stopped at any appropriate stage, forinstance when a predetermined criterion is met. For instance, theiterations may be stopped when two successive solutions A.sup.(n) andA.sup.(n-1) are sufficiently close (the absolute value of δA.sup.(n) isless than a predetermined value). Alternatively, the iterations may bestopped when the reduction in energy represented by the data is lessthan a predetermined amount.

When the source signature has been found to a required precision, it maybe substituted into equation (1). This equation may be truncated at anydesired term so as to remove or reduce the effects of free-surfacemultiples of orders up to and including the highest order of A in thetruncated form of equation (1).

In practice, the updated first order multiple term, D₁.sup.(n), inequation (15), can be kept identical to the initial first order multipleterm, D₁ because primary energy remains unchanged through iterations.

The seismic reflection data required as the input for this technique are"shot gathers", for instance as described with reference to FIG. 1,regularly sampled in space and time. It is assumed that de-ghostingcorrection and subtraction of the direct waves have been performed onthe input data. No assumptions of wavelet structure or phase arerequired in order for the technique to be performed. However, the sourcesignature determined by this technique may be different from the actualsource signature, depending on the effects of preprocessing of the data.

The weighting function W_(A) allows any a priori information about thesource wavelet to be taken into account. For instance, in order toensure that the source spectrum is smooth, this weighting function maybe given as follows: ##EQU4##

The weighting function W_(D) describes errors in the data. It may beused to ensure that certain parts of the data do not dominate thesolution. For instance, this weighting function may be given by thefollowing: ##EQU5## where k_(T) is about one tenth of the maximum wavenumber. Such a function reduces the influence of wave numbers in theneighbourhood k_(g) =0 (i.e. near the zero offset data) on the solution.

The derivation of the present technique has been exemplified for a twodimensional earth model. However, this can readily be modified to athree dimensional earth model with point sources. This simply requiresthat, in equations (9) to (16), the wave numbers k_(s) and k_(g) arereplaced by vector wave numbers ##EQU6##

Further, the inverse source signature can be derived in the (x-w)domain. In this case, equation (8) becomes: ##EQU7## where

    N(ω)=∫dx.sub.g ∫dx.sub.s.D.sub.0 (x.sub.s, x.sub.g, ω)D.sub.1 *(x.sub.s, x.sub.g,ω)               (21)

and

    Q(ω)=∫dx.sub.g ∫dx.sub.s.D(x.sub.s,x.sub.g,ω)W.sub.D (x.sub.s,x.sub.g,ω)D.sub.1 *(x.sub.s,x.sub.g ω)(22)

The fields D₀ (x_(s),x_(g),w) and D₁ (x_(s),x_(g),w) are the temporalFourier transforms of D₀ (x_(s),x_(g),t) and D₀ (x_(s),x_(g),t),respectively. The weighting functions W_(D) (x_(s),x_(g),w) and W_(A)(w,w') must be suitably chosen for this domain.

If, in addition to time variations, the source signature varies withangle or shot position, then the technique described hereinbefore can beperformed for A(w,k_(s)) instead of A(w) in the (w-k) domain or forA(w,x_(s)) instead of A(w) in the (x-w) domain.

FIGS. 2 to 13 illustrate the use of this technique on synthetic datarelating to an acoustic model of the earth comprising two homogeneoushorizontal layers overlying a homogeneous half space. The depth of thetwo horizontal interfaces are 120 and 285 meters. With increasing depth,the velocities of the layers and the half space are 1500 meters persecond, 2000 meters per second, and 2500 meters per second. A syntheticshot point gather was generated comprising 83 receivers with a spacingof 12.5 meters.

FIG. 2 shows the synthetic data as a graphical representation ofrecorded receiver data against receiver offset i.e. in the usual way forrecorded seismic data. The data contain two primaries P1 and P2 butseveral free-surface multiples. This corresponds to the first term D0 inequation (1).

FIG. 3 shows the first two terms D₀ +AD₁ of equation (1). This shows theeffect of removing first order free-surface multiples. Similarly, FIG. 4shows the first three terms D₀ +AD₁ +A² D₂, thus representing removal ofthe first and second order free-surface multiples. FIG. 5 shows thefirst four terms D₀ +AD₁ +A² D₂ +A³ D₃ i.e. with the first, second andthird order multiples removed. FIGS. 6 to 8 show the source signature orwavelet, the amplitude spectrum, and the phase spectrum, respectively,of the actual source signature (solid lines) and the initial estimatedsource signature (broken lines) before iteration to reduce the effectsof truncation as described hereinbefore. Thus, most of the maincharacteristics of the source have been recovered by the firstestimation. In particular, the phase of the estimated source signaturematches well with that of the actual source signature. The result offree-surface multiple removal using this estimated source signature isshown in FIG. 9. Some residual energy remains, for first order multiplesin particular. FIGS. 10 to 13 correspond to FIGS. 8 to 9, respectively,but illustrate the result after five iterations. The match between theactual and estimated source signatures has significantly improved asshown in FIGS. 10 to 12. As shown in FIG. 13, all of the free-surfacemultiple reflections have been successfully removed so that only theprimaries remain.

It is thus possible to provide a technique which permits the signatureof a seismic source to be determined without requiring any assumptions(about the source amplitude or phase or about the subsurface except thevelocity of water in the marine survey) and without requiring anyspecial configurations of hardware. Further, it is possible to provide atechnique based on such source signature estimation to allow the effectsof free-surface multiple reflections to be reduced or removed fromseismic reflection data.

The following situations generate multiple reflected energy in a mannersimilar to that of the marine seismic example described in detailhereinbefore. In such situations, the multiple energy can be used toestimate the source signature which can then be used for the removal orreduction of the multiple energy.

(1) In vertical seismic profiling (VSP) applications, energy sources arelocated on the surface and receivers are located in the subsurface alonga bore hole. Geophones closer to geological horizons can provide a moredetailed stratigraphic picture.

(2) Measurements can be derived from acoustic sources and receivers in abore hole. Acoustic sources in a bore hole may be used to determine (a)the corrosion of the outer surface of the pipe, (b) the third interfacei.e. the bore hole wall, and (c) formation properties near the borehole. These technologies are relevant to well development for drillingto production. Receivers are either coincident with the source or offsetfrom it along the bore hole. The frequency range and radiation patternof the source differ with different applications. In applications (b)and (c), reverberation of the wave in the pipe can have a deleteriousand disadvantageous effect.

(3) In ground penetrating radar, electromagnetic reflection recorded onthe surface of the earth is used to estimate near-surface properties.

(4) In medical imaging and non-destructive material evaluationapplications, a sample is irradiated with energy, such as x-rays orultrasound, in order to determine the internal structure of the samplenon-invasively.

We claim:
 1. A method of estimating a wavelet representing output froman energy source comprising:constructing an inverse scattering serieswhich represents multiple reflected energy and which comprises apolynomial in an inverse of a source signature and has recorded data asthe first term thereof; selecting a sub-series comprising two terms ofthe inverse scattering series, the two terms comprising the first termand an Nth order term where N is a positive integer; and estimating avalue of the inverse signature which substantially minimizes the energyrepresented by the sub-series.
 2. A method of estimating a waveletrepresenting output from a marine seismic energy sourcecomprising:constructing an inverse scattering series which representsmultiple reflected energy and which comprises a polynomial in an inverseof the signature and has recorded marine seismic data as the first termthereof; selecting a sub-series comprising two terms of the inversescattering series, the two terms comprising the first term and the Nthorder term where N is a positive integer equal to 2; and estimating avalue of the inverse signature which substantially minimizes the energyrepresented by the sub-series, comprising analytically solving anexpression for energy minimization of the sub-series to determine theinverse signature.
 3. A method as claimed in claim 2, furthercomprising:replacing the first term of the sub-series by the precedingvalue of the sub-series; replacing the inverse signature in thesub-series by an inverse signature correction to form a modifiedsub-series; determining the value of the inverse signature correctionwhich substantially minimizes the energy represented by the modifiedsub-series; and adding the inverse signature correction to the inversesignature.
 4. A method as claimed in claim 3, in which the steps ofreplacing the first term and replacing the inverse signature arerepeated until a predetermined criterion is achieved.
 5. A method asclaimed in claim 4, in which the predetermined criterion is that theinverse signature correction is less than a predetermined value.
 6. Amethod as claimed in claim 4, in which the predetermined criterion isthat the absolute value of the difference between consecutive values ofthe energy represented by the sub-series is less than a predefinedvalue.
 7. A method of determining a signature of an energy sourcecomprising:constructing an inverse scattering series which representsmultiple reflected energy and which comprises a polynomial in an inverseof the signature and has recorded data as the first term thereof;selecting a sub-series comprising two terms of the inverse scatteringseries, the two terms comprising the first term and an Nth order termwhere N is a positive integer; and determining a value of the inversesignature which substantially minimizes the energy represented by thesub-series and wherein the recorded data comprise transformed recordeddata.
 8. A method as claimed in claim 1, in which the recorded datacomprise recorded data preprocessed to remove or reduce direct wave andghost effects.
 9. A method of determining a signature of an energysource comprising:constructing an inverse scattering series whichrepresents multiple reflected energy and which comprises a polynomial inan inverse of the signature and has recorded data as the first termthereof; selecting a sub-series comprising two terms of the inversescattering series, the two terms comprising the first term and an Nthorder term where N is a positive integer; and determining a value of theinverse signature which substantially minimizes the energy representedby the sub-series, in which the recorded data comprise marine seismicdata.
 10. A method as claimed in claim 1, in which the recorded datacomprise prestack seismic data.
 11. A method of reducing the effects ofmultiple reflected energy in recorded data, comprising:selecting thefirst M terms, where M is an integer greater than 1, of an inversescattering series which represents free-surface multiple reflections andwhich comprises a polynomial in an inverse of an energy sourcesignature; estimating a wavelet representing output from an energysource by a method as claimed in claim 1; and substituting thedetermined signature into the selected first M terms to form data inwhich the effects of free-surface multiple reflections are reduced.